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3.0 SUMMARY 


The preparation of geometric data for input to three-dimensional potential 
flow programs is a very tedious, time-consuming (and therefore expensive) task. 
This report describes a geometry package that automates and simplifies this 
task to a large degree. Input to the computer program for the geometry pack- 
age consists of a very sparse set of coordinate data, often with an order of 
magnitude fewer points than required for the actual potential -flow calculations 
Isolated components, such as wings, fuselages, etc. are paneled automatically, 
using one of several possible element distribution algorithms. Curves of inter 
section between components are calculated, using a hybrid curve-fit/surface-fit 
approach. Finally, intersecting components are repaneled so that adjacent 
elements on either side of the intersection curves line up in a satisfactory 
manner for the potential -flow calculations. The geometry package has been 
incorporated into the NASA Langley version of the 3-D lifting potential -flow 
program and it is possible to run many cases completely (from input, through 
the geometry package, and through the flow calculations) without interruption. 
Use of this geometry package can significantly reduce the time and expense 
involved in making three-dimensional potential -flow cal culations . 



4.0 INTRODUCTION 


With the advent of modern high-speed computers, the aerodynamicist 
gained the ability to study the effects of configuration changes in great 
detail. Whereas before he might have been able to estimate the effects of 
changes in such gross geometric parameters as sweep and aspect ratio on the 
total lift of an isolated wing, he now can accurately calculate the effects 
of changes of quite small details of the shape of very complex configurations 
(wing-fuselage-nacelle-pylon cases, for example) on not only lift but also on 
spanwise and chordwise load distributions, local pressures, flow angles, etc. 
This gain in computational ability has made the preparation of input data for 
describing configuration geometry a very tedious, time-consuming task. General 
three-dimensional potential -flow programs, such as the Hess program, 

(references 1 and 2) require the geometry to be input as the coordinates of 
sets of points, which are grouped to form quadrilateral surface elements. 

Many more elements are generally required to obtain accurate flow solutions 
than the number required for adequate pictorial representation. For example, 
the simple trapezoidal wing of figure 1, which can be represented fairly 
accurately by twenty or so elements, requires on the order of two hundred 
elements to obtain a reasonable flow solution. More complex configurations, 
such as the wing- fuselage case of figure 2, might use nearly a thousand. A 
case which was recently run by personnel of NLR, Amsterdam, using 1780 
elements to represent an external -store configuration is shown in figure 3 
(reference 3). 

There is a clear need to automate as much as possible the preparation of 
the input coordinates, both in order to reduce the number of points input to 
the programs and in order to relax some of the restrictions on how the points 
must be distributed. The main difficulty in doing this is that there are so 
many logically different cases to consider. For example, a wing has different 
spacing requirements than a fuselage; a fuselage has different requirements 
than a nacelle, or a pylon, etc. Moreover, bodies of similar type may require 
different numbers and distributions of points due to the proximity of other 
bodies. If there are intersections between bodies, additional requirements 
are imposed on the point distributions. A geometry package must be very 
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Figure 3. Surface elements used by NLR, Amsterdam, for an external store 
configuration (1780 elements). 

very flexible to apply to so many different cases, but it should not be so 
flexible that it becomes cumbersome to use. 

This report describes a geometry package developed for use with the 
potential -flow program of references 1 and 2 and for extensions of the method 
which may replace it in the future. Since the requirements of future pro- 
grams are mostly speculative at the present, compatibility with the present 
program is emphasized. 

The geometry package has been incorporated into the NASA-Langley version 
of the potential -flow program. This program accepts input either in the 
original input format described in reference 4 (with minor modifications) or 
in the format of the program described in reference 5. With a small amount 
of additional input to control the geometry package, it is possible to run 
a number of frequently occurring cases completely without human intervention. 
For example, isolated wings or bodies may be input to the program with a 
minimum number of points and the geometry package will augment and redistribute 
the points to the number and distribution specified by the user. These 
coordinates can then be punched on cards for use with automatic plotting 
programs to inspect the results before proceeding, or the potential-flow 


program can proceed immediately to analyze the flow. Options are provided to 
allow the user to tailor the distributions to the needs of his particular 
problem, so It should also be possible to run cases having multiple non- 
intersecting wings or bodies with no human intervention. More difficult 
cases, involving intersecting components, can also be treated by the geometry 
package. Some of the simpler cases involving only two components, such as a 
wing/fuselage case, can also be run completely without intervention, but more 
complex cases almost certainly should be checked before proceeding with the 
flow calculation. Some cases cannot be completely handled by the geometry 
package, but in these cases, use of the geometry package to augment the point 
distributions and to calculate intersection curves still results in a 
significant reduction in the effort required to prepare the coordinate data. 
An outline of the major features of the geometry package is given below. 

o Significant reduction in effort required to input a case to the 
potential -flow program. 

o Two separate modes of program input available. 

o Complete compatibility with the potential -flow program — allowing many 
cases to be run completely without interruption between the geometry 
package calculations and the potential -flow calculations. 

o Paneling of isolated components 

o Allows very sparse input coordinate data, 
o Provides output coordinate data suitable for potential -flow 
analysis. 

o Uses independent cubic curve-fits for interpolation in two 
directions on surfaces (N-lines and M-lines). 
o Provides six options for the point distributions on N-lines 
and four options on M-lines. 

o Allows all N-lines on a component to lie in parallel planes, 
if desired. 

o Calculation of intersection curves 

o Requires a distinction between intersecting and intersected 
components. 
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o Intersecting components represented by their M-lines only, 
o Intersected components represented by three-dimensional surface 
f i ts . 

o Intersection curves defined by arrays of intersection points 

between the M-lines of the intersecting components and the fitted 
surfaces of the intersected components. 

o Repaneling of intersecting and intersected components 

o Intersection curve made an N-line of the intersecting component, 
d Other N-lines on intersecting component shifted to restore a 
smooth distribution. 

o Extra strip inside intersected component automatically generated, 
if desired. 

o Intersected components repaneled to insure that elements on 
adjacent components line up along the intersection curves, 
o Simpler repaneling options also provided. 

The remainder of this report documents the theory and operation of the 
geometry package. Section 6 defines the geometric terms used in the later 
sections and discusses the basic philosophy of the geometric input to the 
method. Section 7 describes the paneling of isolated, non-intersecting 
components, including the options available, the applicability and limitations 
of the options, user requirements, methods used, and sample results. Section 
8 covers the method of calculating curves of intersection between components, 
including the theory, restrictions, and verification cases. Finally, 
section 9 describes the methods that have been provided for repaneling 
components after having calculated the curves of intersection between them. 
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5.0 SYMBOLS 


A 

A,B,C,D 

A,b,C,D, 

E.F.G.H, 

I,J,K,L, 

M,N,0,P 

c 

d 


G 

i 

K 

k 

z 

M 

M-line 

N 

N-l ine 

P 

Pi 

S 


Matrix of surface-fit coefficients in algebraic form. 

Coefficients of the equation of a plane. 

Algebraic surface-fit coefficients. 

Local value of the chord of a component. 

Straight-line distance between adjacent points on a curve. 

Also used for the distance from a point to a plane. 

Functions used to express a cubic curve in terms of the functional 
values and first derivatives at the ends of the curve only. 

Dummy variable used to express the general form of a function used 
with several different variables. 

Matrix of surface-fit coefficients in geometric form. 

General subscript used in a variety of ways. 

Parameter used in specifying the distribution of N-lines on a 
component. 

Curvature at a point on a curve. 

Subscript denoting variables associated with a line or curve. 

Matrix of constants used in converting surface-fit coefficients 
from geometric to algebraic form. 

Curve on a component, generally running spanwise on lifting 
components and in the axial direction on nonlifting components. 

Total number of defining points on a curve. 

Curve on a component, generally running chordwise on lifting 
components and in the circumferential direction on nonlifting 
components. 

Total arc length of a curve. 

Normalized point number of the point having index i. 

Arc length along a curve. 
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Sg,S-| Arc lengths at beginning and end of a curve. 

s Subscript denoting a surface. 

T Superscript denoting the transpose of a matrix. 

u,w Parameters used in calculating surface fits. 

x,y,z Coordinates of a point in a Cartesian coordinate system. Also 
subscripts referring to these coordinates. 

3 Angle around a circle circumscribed about an airfoil section, used 

in determining the cosine point spacing distribution. 

e.j Angle of a curve at the defining point having index i. 

e m . Angle of a straight-line segment of a curve beginning at the 

defining point having index i. 
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6.0 NOMENCLATURE AND ARRANGEMENT OF INPUT POINTS 


Before describing the geometry package itself, it is necessary to 
discuss the general scheme for inputting points and for ordering the elements 
that they form and to define some of the terms which are used frequently 
throughout the remainder of the report. Most procedures and definitions are 
identical to those described in references 1 and 2, but some (for example, 
the definition of "ignored elements") have been changed slightly. Reference 
should be made to figure 4 to clarify the discussion which follows. 


STRIP ON 
INTERSECTED 
COMPONENT 


STRIP ON 
INTERSECTING 
COMPONENT 


ELEMENT ON 
INTERSECTING 
COMPONENT 


N-LINES ON 

INTERSECTED 

COMPONENT 



M-LINES ON 

INTERSECTED 

COMPONENT 


ELEMENT ON 
INTERSECTED 
COMPONENT 
N-LINES ON 
INTERSECTING 
COMPONENT 


M-LINES ON 

INTERSECTING 

COMPONENT 


Figure 4. Definition of frequently used terms. 
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A complete configuration (such as a wing-fuselage-nacelle-pylon case) 
is assumed to be constructed of a number of components, each of which is a 
set of associated points. Normally a single component is used to represent 
a complete body (such as a wing or a fuselage), but any number of components 
per body is allowed. There are two types of components — nonlifting and 
lifting. Nonlifting components, such as fuselages or other blunt-ended bodies, 
are represented by source distributions over their surfaces and hence have no 
circulation. Lifting components, such as wings or other bodies with sharp 
trailing edges, are represented by both surface source and dipole distribu- 
tions. Circulation about any section of a lifting component is adjusted in 
such a way as to satisfy the trailing-edge Kutta condition. Lifting 
components also have associated dipole sheets which represent trailing vortex 
wakes. Points on the wakes must be input to the program, as well as points 
on the bodies, and are considered to belong to the same components as the 
associated body points. 

Each component consists of a set of points which can be connected in 
such a way as to form a network of intersecting lines called N-lines and 
M-lines. N-lines on lifting components are the lines running approximately 
in the chordwise direction. They divide a wing, for example, into a number 
of distinct sections. M-lines are the lines connecting corresponding points 
on the N-lines. On nonlifting components N-lines are generally those 
surrounding the major axis of the body (if it is possible to define an axis 
at all). In principle, the roles of the N-lines and M-lines on nonlifting 
components can be reversed without adversely affecting the execution of the 
program, but in most portions of the geometry package, it is assumed that the 
N-lines play the role described above. The potential- flow program does not 
require that all N-lines on nonlifting components have the same number of 
points. However since having variable numbers of points per N-line makes it 
impossible to define M-lines, most applications of the geometry package require 
each N-line to have the same number of points as each other N-line and each 
M-line to have the same number of points as each other M-line. N-lines may 
not cross other N-lines (though they may touch at a point); M-lines also may 
not cross other N-lines (though they may also touch at a point). 
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Points must be ordered so that all points on the first N-line (including 
wake points after the body points for lifting components) are input con- 
secutively, followed by all points on the second N-line and so on. The first 
N-line input may be the one at either extreme of the component, but the 
choice determines the order of input of points along the N-lines. In general, 
the order must be such that the negative of the cross product of the vector 
from one point on an N-line to the next point on the N-line with the vector 
from a point on the N-line to the corresponding point on the next N-line results 
in a vector which is directed to the exterior of the component. This requirement 
may be satisfied on a wing, for example, by ordering the N-lines from tip 
to root and ordering points on the N-lines from the trailing edge along the 
lower surface to the leading edge and back along the upper surface to the 
trailing edge. On a fuselage the requirement is satisfied by arranging the 
N-lines from front to back and the points on each N-line increasing counter- 
clockwise (looking aft). The requirement is also satisfied by reversing both 
the order of the N-lines and the order of the points on the N-lines. However 
parts of the geometry package require that N-lines on fuselages start at the 
front of the body, or more generally, that N-lines on nonlifting components 
start at the end farthest from the first M-line of any other components which 
intersect the component. 

The area between adjacent N-lines on a component is designated a strip. 

Each strip on a lifting component has one characteristic value of the dipole 
distribution and one location where the Kutta condition is satisfied. On 
lifting components, it is possible to specify that a strip have a dipole 
distribution but no source distributions and no boundary conditions on any 
of its elements. Such a strip, called an extra strip, is useful for avoiding 
the abrupt ending of a dipole sheet (which would result in a concentrated 
vortex) along the curve of intersection of two components and for controlling 
the behavior of the dipole sheet near wing-tips when the piecewise linear 
vorticity option is used. Refer to reference 1 for more details. 

The area (generally quadrilateral) between adjacent M-lines on a strip 
is termed an element. Each element has one control point where the boundary 
conditions are satisfied and one characteristic value of the source 
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distribution. It is possible to designate some elements to be ignored elements. 
These elements do not have source distributions and no boundary conditions 
are applied to them. If the component is lifting, ignored elements do have 
dipole distributions, however. References 1 and 2 allow for ignored elements 
to be defined only for lifting components, but this restriction has been 
lifted in the present program. 
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7.0 PANELING OF ISOLATED COMPONENTS 


7.1 General Features of the Paneling Method 

The first operation performed by the geometry package is to panel 
(i.e., distribute the elements) the components as isolated bodies, whether or 
not any of them intersect. Although the resulting element distributions on 
intersecting components may not be useful for analyzing complete configu- 
rations, they serve as a starting point for determining the final distributions 
as well as allowing a configuration "build-up" to be performed (i.e., the 
successive addition of components can be performed, in order to determine 
interference effects). It is assumed that each component is completely 
independent of all others. Since a single body may be composed of more than 
one component, however, and since the close proximity of another body may 
modify the element spacing requirements of a component, it is not always true 
that individual components are completely independent. Sufficient flexibility 
has been designed into the paneling method to allow points to be matched where 
components meet (making their M-lines continuous) and to allow the user to 
specify whatever distribution he deems appropriate to account for the 
proximity of other components. 

The paneling of an individual component is accomplished in two steps. 

First the points on the initial N-lines are augmented in number and redistri- 
buted according to the number and the spacing algorithm specified by the user. 
For this calculation, the N-lines must be roughly chordwise on lifting 
components and, if the later repaneling of intersected components is required, 
the N-lines must be roughly circumferential about an axis on nonlifting 
components. Points in the wake of an N-line of a lifting component are 
distributed independently of the body points. Provision has been made to 
allow the user to input different numbers of points on the initial N-lines 
of a component. The process of redistributing points along the N-lines makes 
the numbers of points equal, thus allowing M-lines to be formed. The second 
step is to augment and redistribute the N-lines according to the number and 
algorithm specified by the user. This is done by augmenting and redistributing 
points along each M-line. 
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The process of redistributing points along either N-lines or M-lines 
requires a method of interpolating along general curves in space. Since 
these curves (particularly the N-lines) are not usually monotonic in either 
x,y, or z coordinates, some other parameter must be used as the independent 
variable of interpolation. In the present method, the parameter chosen is 
the arc length along the polygon formed by connecting straight-lines between 
adjacent points. In the remainder of this report the term "arc length" 
always refers to this straight-line approximation. Interpolation for a point 
on a cutfve requires three separate interpolations, one for each coordinate. 
Each separate interpolation is accomplished in two steps. In the first step, 
numerical differentiation of the dependent variable with respect to the arc- 
length is performed at each defining point on the curve. In the second step, the 
values of the function and its derivatives at the ends of the segments of the 
curve are used to derive the coefficients of cubic interpolating polynomials. 

The numerical differentiation is done using a "weighted-angle" approach 
(see figure 5). In this approach, the lengths (d^) and angles (e m .j) of 
straight-line approximations to segments of the curve are first calculated. 

The angle at the midpoint of each segment is assumed to equal the angle of 
the straight-line approximation to the segment. The angle at any of the 
given points on the curve (e-j) is then determined by taking an average of the 
angles of adjacent segments, weighted by the distances to the midpoints of the 
segments, 


e i = ( d i em-i.-i + e mi )^ d i + d i-l^ ( 7 -l-D 

Angles of the first and last points on a curve are determined by extrapolation, 
assuming constant curvature between the first two and' last two points. The 
derivatives at the given points are then found by taking the tangents of the 
calculated angles. 

In order to facilitate the treatment of components having sections 
bounded by straight lines, such as a wing with a kink in its trailing edge 
or a fuselage with a cylindrical midsection, a multisection component option 
is provided in the numerical differentiation procedure. With this option. 
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Figure 5. Illustration of numerical differentiation procedure. 

the curve is divided into several sections. At the ends of the sections, 
extrapolations are used to determine the derivatives, just as they are at the 
ends of the complete curve. Sections consisting of just two points are 
represented by straight lines. The slopes of these straight lines then 
determine the derivatives at the ends of adjacent sections. 

The values of the function and its derivative at each end of a curve 
segment constitute four pieces of information which can be used to determine 
the coefficients of a cubic polynomial approximating the curve. The form of 
the polynomial is 

f(S) = f (S 0 ) + (S - S 0 ) f ' (S 0 ) + }(S - S 0 ) 2 f"(S 0 ) + -g-{S - S 0 ) 3 f M, ($ 0 ) 

(7.1.2) 

where S is the independent variable (arc length ranging from values s 0 
to S-|, f represents the dependent variable (x, y, or z coordinates), and 
primes denote differentiation. The form of its derivative is 

f'(S) = f'(S 0 ) + (S - S 0 ) f"(S 0 ) + \ (S - S 0 ) 2 f'feo) (7.1.3) 
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Given the values of f and f' at S Q and Sj, the simultaneous solution 
of equations (7.1.2) and (7.1.3) yields the values of f"(S 0 ) and f'"(S o ). 
Equation (7.1.2) can then be used to determine the value of the function f(S) 
at any value of S within the given curve segment. 

This curve-fit method does not insure continuity of the second deriv- 
ative of the function and thus is not a cubic spline fit in the usual sense 
(reference 6). It was chosen rather than a spline method because of its 
consistently superior results in several test cases. For example, figure 6 
shows a comparison between this method and a true cubic spline method. 

The interpolated coordinates found by the present method are considerably 
less wavy than those calculated using the spline method. Both methods 
fail to represent the shape accurately in the after region of the body, 
because of the higher local curvature and the proximity of the ends of the 
curve. Other comparisons involving airfoils and more general shapes also 
showed smoother results for the present method than the spline method, in 


INPUT POINTS 

INTERPOLATED POINTS (PRESENT METHOD) 
INTERPOLATED POINTS (SPLINE METHOD) 



Figure 6. Comparison of curve-fit methods. 
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spite of the alleged proof in reference 6 that, in general, a spline produces 
the smoothest of all possible curve fits. The contradiction of this proof 
may possibly be due to the unusual character of the independent variable (it 
being defined as a quantity which varies smoothly along a nonsmooth curve) 
or it may be due to the presentation of the results by a graph of one 
interpolated result (the z-cocrdi nates) as a function of other interpolated 
results (the x-coordinates). 

7.2 Distribution of Points Along N-Lines 

Each of the options for distributing points along N-lines (except the 
trivial option of leaving the initial distribution unchanged) requires an 
array of normalized arc lengths which applies to every N-line of the com- 
ponent under consideration. The formulas for these arc-length distributions 
are given below. Given the specified distributions, the method calculates 
the distribution on each initial N-line and interpolates each coordinate 
independently to determine the values at the desired locations. In some 
cases all points on an N-line coincide (as at the end of a pointed body, for 
example). Then the method does not attempt to interpolate, but simply 
provides the specified number of points to that N-line. The following options 
are available for distributing points on N-lines: 

1. Input distribution, unaltered 

2. Input distribution, augmented in number 

3. Constant increments in arc length 

4. Constant increments on the superscribed circle (cosine spacing) 

5. Curvature-dependent distribution 

6. User-specified distribution 

7.2.1 Input Distribution, Unaltered 

With this option the method does no interpolation. It should be used 
whenever the initial distribution already contains a sufficient number of 
properly spaced points. If the N-lines of the component under consideration 
are to be redistributed, or if the component is involved in an intersection with 
another component, then the number of points input must be the same on each 
component of the N-line. In addition, the distribution of points on adjacent 
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N-lines should be fairly similar so that when corresponding points on N-lines 
are connected, the resulting M-lines make smooth curves. This is the only 
option which may result in different numbers of points and different 
distributions of points among the N-lines of a component. Although no 
interpolations are required, some calculations are still necessary to deter- 
mine the derivatives with respect to arc length of the coordinates at each 
point on each N-line. These are needed in later calculations to determine 
intersection curves. 

7.2.2 Input Distribution, Augmented in Number 

This option results in a distribution which is similar to the initial 
distribution, but contains a different number of points. The initial distri- 
bution used in the calculation is the distribution for the first N-line on 
the component, unless all points coincide on the first N-line. In that case 
the distribution on the second N-line is used. The method works by defining 
a "normalized point number", p^ , ranging from zero to one. 

Pi = (1 - 1 )/ (N - 1) (7.2.1) 

where i is the index of the point and N is the total number of points. 
Arrays of normalized point number are formed for both input and output dis- 
tributions. The output arc-length distribution is determined by interpolating 
the curve of input arc length versus input normalized point number to the 
output values of the normalized point number. This procedure is illustrated 
in figure 7. Suitable cases for this option include both lifting and non- 
lifting type components. It is especially useful for precisely controlling 
the desired distribution without having to load a large number of points and 
without having to determine in advance what the numerical values of the 
arc lengths are. Typical results for a section of a supercritical wing are 
shown in figure 8. 
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NORMALIZED POINT NUMBER - [(i — 1 )/(N — 1 )] 

Figure 7. Method of distributing points on an N-line — input distribution 
augmented. 


+ INPUT POINTS (N = 15) O OUTPUT POINTS (N = 25) 



Figure 8. Point distribution on a supercritical wing section — input distribu- 
tion, augmented. 
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7.2.3 Constant Increments in Arc Length 

This option results in points distributed uniformly according to the 
formula 

Si/P = (1 - 1)/(N - 1) (7.2.2) 

where S is the arc length measured from the first point and P is the 
total arc length around the perimeter of the N-line. It should be used for 
smooth bodies which do not have large variations in curvature (for example, 
the cylindrical fuselage section illustrated in figure 9. It should usually 
not be used for wing sections, since, for a reasonable number of points, there 
would be too few points to accurately define the shape near the leading edge 
or other high-curvature regions and too coarse spacing at the trailing edge 
for the Kutta condition to be accurately applied. 



Figure 9. Point distribution on a cylindrical fuselage section — constant 
increments in arc length. 

7.2.4 Cosine Spacing 

Cosine spacing is a commonly-used distribution, dating from some of the 
classical two-dimensional theories (reference 7). Its name derives in an 
obvious way from the usual formula for the spacing 
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(7.2.3) 


x-j/e = \ (1 + cos[(i - l)n/(N - 1)]} 

where c is the chord, which is assumed to be parallel to the x-axis. The 
significance of this method, as illustrated in figure 10, is that points are 
chosen having the x-coordinates corresponding to equal increments in angle 
around a circle circumscribed about the N-line with leading and trailing edges 
touching the circle. This results in very fine spacing near leading and 
trailing edges and coarse spacing in regions away from these edges. It is a 
useful distribution option for wing sections but would probably not be 
commonly used for nonlifting components. 



Figure 10. Method of distributing points on an N-line — cosine distribution. 

Equation (7.2.3) is not directly useful for determining the specified 
arc-length distributions which the method of this geometry package requires. 

To do this, the wing section is first scaled to unit chord, translated to 
put the origin of the coordinate system at the leading edge, and rotated to 
make the chord parallel to the x-axis. The angle (e) about the superscribed 
circle is then calculated for each point on the N-line by solving the equation 

x i/c = (1 + cos 3) (7.2.4) 
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using the transformed values of the x-coordinates. The arc-length distribution 
corresponding to the input points is also calculated. For both the angle and 
arc-length calculations, only the first N-line on a component is used (unless 
its length is zero, in which case the second one is used). Given the values 
of the arc-length distribution and the angle distribution for the input 
coordinates, and the desired values of the angle distribution for the output 
coordinates (uniform increments), the arc-length distribution for the output 
coordinates is determined by interpolation. 

Results of this option for the same supercritical wing section used 
previously are shown in figure 11. 



Figure 11. Point distribution on a supercritical wing section — cosine 
distribution. 

7.2.5 Curvature — Dependent Distribution 

The accuracy and computational efficiency of the potential -flow method 
are both improved by using fine spacing in regions of high curvature and 
coarser spacing in regions of low curvature. This option produces a point 
distribution which depends on the curvature distribution along the first 
N-line of the component (or the second if the length of the first is zero). 

If the component is of lifting type, the point distribution is also made to 
be a function of the proximity to the trailing and leading edges in order to 
allow the trailing-edge Kutta condition to be accurately satisfied and to 
keep the variation of the spacing smooth in the leading-edge region where the 
curvature varies rapidly. 
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The following relationship between the curvature and the spacing increments 
is employed: 


4s i - - A W As ra ax )(1 -VW + 4s n,in /4 W 4s n,ax 


where as. 


(7.2.5) 

is the ith increment in normalized arc length between adjacent 


points, k.j is the absolute value of the curvature at the center of the ith 

segment, k, is the maximum absolute value of the curvature on the N-line, 
max 

and 


and as 


min 


AS 


max 


are, respectively, the minimum and maximum allowable 

is 


increments in arc length. In this application, the ratio As m . /as„„„ 

3 rr min max 

specified to be 0.25. 


To implement this relationship requires an iterative procedure. First 

the arc-length and curvature distributions are calculated. Then, if the 

component is of lifting type, the variable which plays the role of the 

curvature in equation (7.2.5) is modified to make its value at the trailing 

edge equal to its leading-edge value and to make it vary linearly from the 

leading and trailing-edge values to the values at 7.5 percent of the N-line's 

perimeter on either side of the leading and trailing-edge points. This 

modification insures that points will be closely spaced near the trailing 

edge and that the spacing will not vary abruptly near the leading and trailing 

edges. To start the iterations, estimates of the values of As max and all 

the values of as^ are required. These are all taken to be equal to the total 

length of the curve divided by the number of segments it contains. Given the 

estimated values of As ma „ and as^, the curvature at the center of each 

segment (modified as described above) is determined by interpolation of the 

curvature versus arc-length relationship and used in equation (7.2.5) to 

update the estimated values of as^ . These values are then scaled to make 

their sum equal to the total length of the N-line and searched to determine 

the value of as . These updated values of as. and as„,,„ are then used 
max i max 

in the next iteration and the process is repeated a maximum of five times or 
until the required scale factor for the lengths of the increments is 
sufficiently close to unity. 

Results of this option for the supercritical wing section used previously 
are shown in figure 12. Results for a section of a nonlifting component are 
shown in figure 13. 
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Figure 12. Point distribution on a supercritical wing section-curvature- 
dependent distribution. 



Figure 13. Point distribution on a section of a nonlifting component— 
curvature-dependent di stri buti on . 

7.2.6 User-Specified Distribution 

For N-lines of unusual shape or when components are located very close 
together, it is possible that none of the options described above will produce 
an adequate point distribution. To account for this possibility, provision 
has been made for the user to specify (as input data) values of the normalized 
arc lengths, which are used directly to determine the augmented coordinates. 

7.3 Distribution of N-lines 

After the points on the initial N-lines have been redistributed, each 
N-line on a component has the same number of points distributed in a similar 
manner. Connecting corresponding points on all N-lines generates a set of lines 
designated M-lines. Since the distribution of points is the same on each N-line, 
the M-lines are smooth and have fairly small curvature. The process of 
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augmenting and redistributing the N-lines is accomplished by augmenting and 
redistributing points along these M-lines. For these calculations, it is 
required that each M-line have a total length greater than zero. 

There are four options available for the distribution of M-lines. They are: 

1. Input distribution, unaltered 

2. Input distribution, augmented in number 

3. Constant increments 

4. User-specified distribution 

The descriptions of all four options are completely analogous to the descrip- 
tions of the corresponding options for distribution of points on N-lines 
described in section 7.2. 

There are two modes of operation of this portion of the geometry package. 

The first (and most common) mode is designated the planar-section mode and the 
second is termed the arc-length mode in the discussions which follow. 

7.3.1 The Planar-Section Mode 

The usual way of inputting data to the potential -flow program is to make 
the N-lines on a component planar cuts perpendicular to some axis of the 
component. For nonlifting components, this convention is just a convenience 
to the user. For lifting components, however, the elements are approximated 
by trapezoids, so that if the N-lines are not parallel, the accuracy of repre- 
sentation of the component by its elements is not as good. 

The planar section mode of operation redistributes the points on the 
M-lines in such a way that all N-lines, except perhaps the first and last ones 
on a component, are parallel. The first and last N-lines generally define the 
planform view of the component, which is not required to be composed only of 
straight lines. The method does not alter these N-lines during this portion of 
the calculations. The distributions described above generally refer to the 
distances between the parallel planes, normalized by the distances between the 
first and last N-lines. At the edges of the component, when the first and last 
N-lines are not planar, the distributions refer to distances between the first 
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point on the N-line under consideration and the first point on the first N-line 
of the component, resolved in the direction perpendicular to the parallel planes 
and normalized by the distance between the first points on the first and last 
N-lines of the component. The orientation of the planes is defined by specify- 
ing (as input data) the direction cosines of a vector perpendicular to the planes 
For this calculation, the sense of the normal vector is not important. However, 
if intersection curves are calculated and the components are subsequently repan- 
eled, then the vector must point in the direction of increasing N-lines. If the 
direction cosines are not specified and the planar section mode is used, default 
values of (1.0, 0.0, 0.0) for nonlifting components and (0.0, -1.0, 0.0) for lift 
ing components are assumed. When the default values are used, N-lines must be 
input from front to back on nonlifting components and from tip to root of lifting 
components. 

The first step in distributing the points on the M-lines is to determine 
the equation of the planes representing the N-lines and the interpolating poly- 
nomials representing the M-line segments. The equations of the planes of the 
N-lines are of form 

Ax + By + Cz + D i = 0 (7.3.1) 

where A, B, and C are equal to the x, y, and z direction cosines of 
the vector normal to the planes (and hence are the same for all planes on the 
component). The values of D. (different for each plane) are given by 

D i = D 1 + K i (D 2 “ D i> (7.3.2) 

where D-j and D 2 are the coefficients of the planes passing through the 
first and last points, respectively, on the first M-line of the component and 
fCj is the value of the distribution parameter for the N-line under consider- 
ation. D-, is determined from the relationship 

D 1 = -(Ax-, + By-, + Cz., ) (7.3.3) 

where x-, , y., , and z., are the coordinates of the first point on the first 
M-line. D 2 is determined in a similar manner, using the last point on the 
first M-line. The interpolating polynomials for the M-line segments are 
determined by the curve-fit procedure described in section 7.1. 
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The second step in distributing points on the M-lines is to determine which 
M-line segments intersect which planes. For this step, it is assumed that M-line 
curvature is small enough so that straight-line approximations to the M-line seg- 
ments may be used. 

Then, knowing which M-line segments intersect which planes and knowing the 
equations of the M-line segments and the planes, the coordinates of each inter- 
section point can he computed. To do this for any point, substitute the expres- 
sions for x,y,z as a function of arc length along the M-line segment (each 
having the form of equation (7.1.2)) into the equation of the plane (7.3.1). This 
gives a single cubic equation with S-j , the arc length along the M-line segment 
to the intersection point, as the only unknown. This is solved by Newton's 
method. Using the intersection point of the plane and the straight-line approx- 
imation to the M-line segment to determine the starting point, the method typi- 
cally requires only three or four iterations to converge to an adequate solution. 

Figures 14 and 15 show the top views of two typical wings both before 
and after redistributing the N-lines. Since the wing in figure 14 is trape- 
zoidal in planform and has a linear twist distribution, only the N-lines at 
the tip and root are needed in the initial geometric representation. Since 
the initial representations of the tip and root sections are already planar 
and all M-lines are straight lines, the planar-section mode and the arc-length 
mode (described in the next section) produce identical results. A more gen- 
eral case, for which the two modes of operation should give very different 
results is shown in figure 15. In this case, the M-lines on the entire out- 
board half of the span are straight lines, but they are curved on the inboard 
portion. Therefore, the multi section component option described in section 7.1 
is used, allowing the linear portion of the geometry to be represented using 
only two N-lines. Figure 16 shows the top view of a fuselage before and after 
redistributing the N-lines. This fuselage has a cylindrical midsection with 
semi -ellipsoidal sections fore and aft. Again, the multisection component 
option allows the linear portion to be represented using only two N-lines. In 
all three of these cases, the redistributed N-lines are equally spaced. 
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INITIAL ELEMENT 
DISTRIBUTION 


ELEMENT DISTRIBUTION 
AFTER REPANELING 



Figure 14. Redistribution of elements on a trapezoidal wing (cosine spacing 
chordwise, constant increments spanwise, planar. section mode). 
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Figure 15. Redistribution of elements on a supercritical wing (cosine spacing 
chordwise, constant increments spanwise, planar-section mode, 
multisection component option). 
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INITIAL ELEMENT DISTRIBUTION 




ELEMENT DISTRIBUTION AFTER REPANELING 



Figure 16. Redistribution of elements on a fuselage (constant increments around circumference and in 
axial direction, planar-section mode, multisection component option). 


7.3.2 The Arc-Length Mode 

In some cases the N-line at the edge of a component may be so highly nonplanar 
that a strip bounded by this N-line and an adjacent planar N-line would vary 
strongly in width, resulting in areas too sparsely covered with elements. An 
example of this is the under-wing pylon shown in figure 17(a). In this figure 
the thickness of the wing has been exaggerated to help illustrate the point. In 
such cases the arc-length mode of distributing N-lines on the component should be 
used. In this mode the specified distribution parameters refer to fractions of 
the total arc lengths of the M-lines, rather than to normalized distances between 
planes. The points on all M-lines of a component are distributed in a similar 
manner, so that the N-lines in the vicinity of the nonplanar edge of the compo- 
nent curve smoothly around to create a more uniform distribution of the elements, 
as illustrated in figure 17(b). In this mode of operation, the method of redis- 
tributing points along M-lines is completely analogous to the method of redis- 
tributing points along N-lines. 



(a) Planar-Section Mode. (b) Arc-Length Mode. 

Figure 17. Comparison of planar-section and arc-length modes of distribution of 
N-lines — strut on a thick wing. 
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8.0 CALCULATION OF INTERSECTION CURVES 


8.1 General Features of the Method 

The second major operation performed by the geometry package is the calcu- 
lation of the curves of intersection between body components. For this calcu- 
lation it is assumed that the M-lines of one of the components (designated the 
intersecting component) pierce the surface of the other component (designated 
the intersected component). The final solution is a set of intersection points, 
one for each augmented M-line on the intersecting component. In general, the 
method calculates the intersection point of a curve and a surface and hence 
requires a surface-fit method as well as a curve-fit method. The surface-fit 
method uses the theory of parametric cubic surface patches originated by Coons 
(reference 8) and extensively developed and applied by a number of other 
investigators (reference 9, for example). The curve-fit method is the same 
one used in the paneling of isolated components and described in section 7.1. 

8.2 Restrictions and Limitations of the Method 

The intersection method is designed to handle cases which occur frequently 
in aircraft applications, such as wing-fuselage intersections, wing-pylon 
intersections, etc. It is not designed for complex cases in which the inter- 
section curve is discontinuous or for cases involving fillets or smooth transi- 
tions from one component to the next. Use of the method is limited to cases 
in which the following restrictions apply: 

1. A distinction can be made between intersecting and intersected 
components . 

2. A component can intersect only one other component and can be inter- 
sected by only one other component. If a body intersects or is 
intersected by several bodies, it must be divided up into several 
components. A single component can intersect one component and be 
intersected by another component, however. 

3. The M-lines on the intersecting component pierce the surface of the 
intersected component at a sharp angle. Surfaces are not tangent 
where they meet. 


32 



4. The M-lines on the intersecting component extend sufficiently far 
into the interior of the intersected component to allow a planar 
representation of the elements on the intersected component to be 
used in the process of searching for the elements which are inter- 
sected by the M-lines. 

5. The intersecting component has at least one N-line which lies entirely 
in the interior of the intersected component. 

6. No M-line on the intersecting component intersects the intersected 
component more than once. 

8.3 Details of the Method of Solution 

The method of calculation of the intersection curve between two components 
requires that one component be assigned the role of the intersecting component 
and the other component be assigned the role of the intersected component. 

This is done by the user in the input data to the geometry package. The inter- 
section curve is then defined by the set of intersection points between the 
M-lines on the intersecting component and the elements on the intersected 
component. In order to define this curve at a sufficient number of points, 
every M-line on the intersecting component (after redistributing and augmenting 
the points on the N-lines) is used. On the intersected component, in principle 
either the original input elements or the elements after augmenting points on 
just the N-lines, or the elements after augmenting points on both N-lines and 
M-lines could be used for defining surface patches which are intersected by the 
M-lines of the intersecting component. Since parametric cubic surface patches 
are used, however, it is not necessary to use the large number of elements that 
exist after augmenting points on both N-lines and M-lines. Since the input 
points may vary in number from one N-line to the next, and since the distribu- 
tions may also be dissimilar, use of the original input elements could result 
in surface patches having much more extreme curvature (and hence less accuracy) 
than use of either of the other two sets of elements. Therefore, the elements 
after augmenting points on the N-lines but before augmenting points on the 
M-lines are used for the formation of surface patches. 

One possible method for calculating intersections of an M-line and a 
surface might operate as follows: Derive mathematical formulas to represent 
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each segment of the M-line and each element of the surface and find the roots of 
the equations obtained by equating coordinate values on the segments and the ele- 
ments. If a particular equation has no real roots, or if it does have real roots 
but the points they represent fall outside the bounds of the element or the seg- 
ment, then the element and the segment have no intersection points. If the equa- 
tion does have real roots and the points they represent do fall within the bounds 
of both the element and the surface, then the roots are the intersection points. 
When all such equations have been solved, the intersection points are grouped, 
somehow, to determine the intersection curves. Such a method would have to solve 
a large number of equations (a number equal to the number of M-line segments 
multiplied by the number of surface elements) and many of these equations would 
fail to have solutions. It would, therefore, be very uneconomical to use. 

A more economical approach, adopted in the present method, first compares 
maximum and minimum coordinate values on the M-line segments and the surface 
elements and eliminates from consideration those combinations which obviously 
contain no intersection points. Then an approximate intersection method, 
with planar element representations, is used to determine which M-line segments 
and surface elements contain intersection points and to get approximate values 
of the intersection point coordinates. Then, mathematical representations are 
obtained only for those elements containing intersection points, and these are 
used, together with the representations of the M-line segments, to calculate 
the intersection points more accurately. Details of this approach are given 
below. 


8.3.1 The Initial Search for M-Line Segments and Surface Elements Containing 
Intersection Points 

The basic information needed to conduct the initial stages of the search 
includes the maximum and minimum values of each coordinate on each element of 
the intersected component, on each strip of the intersected component, on each 
segment of the intersecting M-line, and on the entire intersecting M-line. 
Because it is assumed that the intersecting M-lines extend some distance into 
the body and a planar element representation can be used, it is legitimate to 
assume that the maximum and minimum values on an element must occur at one 
of its corners. A quick comparison of the coordinates of the corner points 
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thus determines the extrema. Extreme values on a strip are determined by com- 
paring the extreme values of each of its elements. Values on the M-line segments 
are determined more accurately, assuming cubic relationships between the coord- 
inate values and the straight-line arc lengths along the M-line (equations (7.1.2) 
and (7.1.3)). Equation equation (7.1.2) to zero (separately for x,y, and z) and 
solving gives two values of arc-length along the M-line (for each coordinate). 

The real solutions which fall within the bounds of the segment are used to calcu- 
late coordinate values which are compared, along with the values at the ends of 
the segment, to' determine the extreme values. Maximum and minimum values on the 
entire intersecting M-line are determined by comparing the extreme values of each 
of its segments. 

Having determined the extreme values of the coordinates on strips and 
elements of the intersected component and on segments and the entirety of the 
intersecting M-line, the searching procedure begins by comparing the extreme 
coordinate values on the first strip of the intersected component with the 
values for the entire M-line. If the range of any of the coordinate values 
on the strip fails to include at least part of the range of the values of the 
same coordinate on the M-line, then no intersection can possibly occur on that 
strip. Successive strips are compared with the M-line until one is found 
which may possibly contain the intersection point. Then each element on the 
strip is compared with the M-line until a possible intersected element is 
found. Then each segment on the intersecting M-line is compared with the ele- 
ment until a possible intersecting segment is found. 

8.3.2 The Final Search for M-Line Segments and Surface Elements Containing 
Intersection Points and the Approximate Determination of the Inter- 
section Points 

Having found an element on the intersected component and a segment of 
the intersecting M-line which share, at least partially, the ranges of their 
x, y, and z-coordi nates , a more careful check on whether or not they inter- 
sect is made. If the element is not already triangular, it is divided into 
two triangles by drawing one of the diagonals across the quadrilateral region. 

Each element is then represented by two planar, triangular subelements. The 
planar coefficients of each subelement are determined, starting from the 
general form of the equation of a plane 
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Ax + By + Cz + D = 0 


(8.3.1) 


by first dividing through by D (illustrating that there are only three 
independent unknowns in equation (8.3.1)) and rearranging to get 

A/Dx + B/Dy + C/Dz = -1 (8.3.2) 

and then substituting the coordinates of the corner points of the subelement 
into equation (8.3.2) to obtain a third-order system of linear equations which 
is easily solved for A/D, B/D, and C/D. A, B, and C are then found by 
multiplying by any arbitrarily chosen (finite, nonzero) value of D. 

A check is then made to see if the ends of the M-line segment lie on 
opposite sides of the plane. This is done by computing the directed distance 
from each point to the plane, using the formula 

Ax-, + By-, + Cz, + D 

d = (8.3.3) 

Va 2 + B 2 + C 2 

where d is the distance and (x-j , y-j , z-j ) are the coordinates of the point. 

If the value of d has the same sign for both end points, the segment does 
not intersect the plane of the triangular subelement. If this is true for 
both subelements on the element, the element and the M-line segment do not 
intersect and the next segment is considered. 

When an M-line segment which does cross the plane of one of the subelements 

is found, the point of intersection of the segment and the plane is found using 

the method described in section 7.3.1. For this calculation, the M-line seg- 
ment is again represented by a cubic polynomial. Given the intersection point, 
it is next necessary to check whether or not it falls within the triangular 
region of the subelement. To do this, the present method checks that each side 
of the triangle lies on the same side of both the intersection point and the 
opposite apex of the triangle. To determine that two points, A and B, lie 
on the same side of a line from points 1 to 2, the method takes the cross 
product of the vectors from A to 1 and from A to 2 and the cross 

product of the vectors from B to 1 and from B to 2 and checks whether 

or not the dot product of the two resulting vectors is positive. 
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If the point does not fall within the bounds of either triangular subele- 
ment, the next segment of the M-line is checked. If none of the segments 
intersects the element, the next element is checked in the above manner. If 
no element on the current strip contains the intersection point, the entire 
procedure, including the initial search described in section 8.3.1, is 
repeated for the next strip, and so on, until either an element of the inter- 
sected body and a segment of the intersecting M-line which contain the inter- 
section point are found, or all possible combinations of elements and M-line 
segments are exhausted. If an approximate intersection point is found in 
this manner, it is used as a starting point for the iterative method des- 
cribed below, to determine a more precise intersection point. If the method 
fails to find an approximate intersection point, calculations continue, 
starting with the next M-line on the intersecting component. In this case, 
execution of the computer program is terminated after the intersection calcu- 
lations since the final repaneling method requires the intersection curve to 
be uninterrupted and completely defined. 

8.3.3 Derivation of a Mathematical Representation of the Surface of an 
Element 

In order to compute a more accurate intersection point, it is necessary 
to obtain an equation for the surface of the element. This is done using the 
theory of parametric cubic surfaces described in references 8 and 9. The 
following discussion is a summary of relevant portions of the theory. 

Required geometric quantities at each corner point of the element under 
consideration include the following (as well as the corresponding y and z 
val ues) : 

x ’ V V x uw’ x wu 

where the subscripts indicate differentiation and u and w are the param- 
eters upon which the surface fit is based. On the boundary curves of an 
element, u represents a fraction of the straight-line arc length between 
the two M-line boundary curves and w represents a fraction of the straight- 
line arc length between the two N-line boundary curves. On the interior of 
an element, u and w represent quantities analogous to those defined above 
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but, for elements with compound curvature, the physical interpretation of the 
parameters is less obvious. x u at a corner indicates the derivative of x 
with respect to straight-line arc length along the N-line passing through the 
corner, normalized by the straight-line distance between adjacent points 
on the N-line. x indicates the derivative with respect to straight-line 

W 

arc length along the M-line passing through the corner, normalized by the 

straight-line distance between adjacent points on the M-line. x and x 

uw wu 

indicate cross derivative terms in the directions of the M-line and N-line. 
The first derivative terms (not normalized) were previously required for the 
calculations to redistribute and augment points along the M-lines and N-lines. 
They need only be normalized by the proper element side length to be applic- 
able here. The cross derivative terms x.,,, and x are obtained by numer- 

uw wu 

ically differentiating the unnormalized first derivative terms 
along the M-lines and N-lines, respectively, and then normalizing by the 
product of the lengths of the adjacent element sides. The two cross-derivative 
values at each corner are then averaged, since the equation used to represent 
the surface implies that they are equal. 


These geometric quantities, for x, y, and z and for each corner of the 
element, constitute the so-called "parametric cubic patch coefficients in 
geometric form." The coefficients can be arranged in matrix form as follows 
(shown only for the x coordinates): 


x 00 

X 01 

o 

o 

** 

Xw oi 

X 

o 

x n 

X W 

w 10 

Xw n 

x u 

u 00 


Xuw oo 

Hi 

x u 

10 


x uw 10 

Hi 


(8.3.5) 


In this matrix the u and w subscripts again represent derivatives as above. 
The 0 and 1 subscripts together indicate the corner of the element being 
considered, the first indicating the N-line and the second indicating the 
M-line. The 0 indicates the lower-numbered line on the element and the 1 
indicates the higher-numbered line. Similar matrices are also formed for the 
y and z-coordinate terms. 
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The matrix of coefficients in geometric form, [G x ], and the corresponding 
y and z matrices contain all the information needed to derive algebraic 
expressions for the coordinates at any point on the surface of the element in 
terms of the parametric variables. These expressions are of the following form 
(again shown only for the x coordinates): 


x(u,w) = w 3 [A x u 3 + B x u 2 + C x u + D x ] 

+ w 2 [E x u 3 + F x u 2 + G x u + H x ] 

3 2 . <8-3.6) 

+ w [I x u J + J x u^ + K x u + L x ] 

+ [M x u 3 + N x u 2 + 0 x u + P x ] 

The coefficients of equation (8.3.6) can be grouped to form the so-called 
"matrix of coefficients in algebraic form," [A ] as follows: 

A 


(8.3.7) 


To show how the matrix of coefficients in geometric form may be converted 
into the matrix of coefficients in algebraic form, it is first necessary to study 
the properties of equation (8.3.6) further. On the boundaries of the element 
one of the parametric variables is constant, either zero or one, and the other 
varies from zero to one. For definiteness, assume that the variable w is 
equal to zero. Equation (8.3.6) then reduces to a cubic equation with u as 
the only independent variable. On the opposite side of the element, w is 
equal to one and equation (8.3.6) reduces to another cubic equation, again 
with u as the only independent variable. The entire surface of the element 
can be considered to be a collection of cubic curves with w constant and u 
variable. Each of these curves has the form of equation (7.1.2), which can be 
converted to the form 

f(u) = F,(u)f(0) + F 2 (u)f(l) + F 3 (u)f'(0) + F 4 (u)f'(l) (8.3.8) 


[A x ] = 


A x B x C x 
E x F x G x 
>x J x K x 
M x N x K 
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by solving for the second and third derivative terms in (7.1.2) as described in 
section 7.1. When this is done, the terms F-j(u), F 2 (u), Fg(u), and F^(u) are 
given by 


F n (u) = 2u 3 - 3u 2 + 1 
F 2 (u) = -2u 3 + 3u 2 
F 3 (u) = u 3 - 2u 2 + u 
F 4 (u) = u 3 - u 2 


(8.3.9) 


or in matrix notation 


[F] = CF,(u) F 2 (u) F 3 (u) F 4 (u)] = [u 3 u 2 u 1] • [M] 


where 


[M] 


+2 

-2 

+1 

+1 

-3 

+3 

-2 

-1 

0 

0 

+1 

0 

+1 

0 

0 

0 


In matrix notation, equation (8.3.8) becomes 


(8.3.10) 


(8.3.11) 


F (u ) = [F] ■ [f(0) f(l) f'(0) f'(l)] T (8.3.12) 

= [u 3 u 2 u 1] • [M] • [f(0) f(l) f 1 (0) f'(l)] T 


where T indicates that the transpose of the matrix is to be taken. Now, using 


(8.3.12) 

, the equations ■ 

for 

x(u) 

and 

yuj 

on 

the boundary 

and w = 

1 can be written as 






x(u,0) = 

£u 3 

u 2 

u 1] • 

CM] ■ 

■ »00 

x 10 

\o 

x u ]T 
U 10 


x(u,l) = 

[U 3 

u 2 

U 1] • 

[M] ■ 

' [x 01 

X 11 

x “ 01 

x u ]T 
U 11 


x w (u»0) = 

[U 3 

u 2 

u 1] • 

CM] ■ 

' [Xw oo 

X w 

10 Xuw 00 Xuw 10 


x w (u,l ) = 

[U 3 

u 2 

U 1] • 

CM] - 

’ [ \l 

x w, 

n Hi Hi 


(8.3.13) 
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or more compactly 


[x(u,0) x(u,l) x w (u,0) x w (u,l)] = [u 3 u 2 u 1] • [M] • [G x ] (8.3.14) 

Equation (8.3.12) can also be used to derive the expression 

x(u,w) = [w 3 w 2 w 1] • [M] • [x(u,0) x(u,l) x w (u,0) x w (u,l)] T (8.3.15) 

or equivalently 

x(u,w) = [x(u,0) x(u,l) x w (u,0) x w (u,l)] • [M] T • [w 3 w 2 w 1] T (8.3.16) 

Combining (8.3.14) and (8.3.16) gives 

x(u,w) = [u 3 u 2 u 1] • [M] • [G x ] * [M] T • [w 3 w 2 w 1] T (8.3.17) 

Since equation (8.3.6) can be rearranged and written as 

x (u ,w ) = [u 3 U 2 u ]] * [A x ] * [w 3 w 2 w 1 ] T (8.3.18) 

the matrix of coefficients in geometric form is converted to the matrix of 
coefficients in algebraic form by the operation 

[A x ] - [M] • [G x ] • [M] T (8.3.19) 

Expressions for the y and z coordinates as functions of the parameters u 

and w are obtained in a similar manner, using the matrices of geometric 

quantities [G ] and [G ]. 

y ^ 

Including coefficients for all three coordinates, there are forty-eight 
parametric cubic quantities associated with each element. Except for the 
coordinates themselves (twelve of the forty-eight quantities), these quanti- 
ties are generally not shared by adjacent elements. However, these coeffici- 
ents are all derived from geometric data which is continuous from element to 
element and, therefore, is shared by adjoining elements. Therefore, a large 
reduction in the computer program storage requirements can be made, at a 
small expense of additional computation time, by storing all data in the 
(unnormalized) geometric form and converting to the algebraic form of equation 
(8.3.6) only when actually needed. Therefore, the coordinates, their deriva- 
tives along the N-lines, and M-lines, their cross -derivatives, and the arc 
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lengths between adjacent points are the only geometric data used in the present 
method. 

8.3.4 Computation of More Precise Values of the Coordinates of the Inter- 
section Point 

At this point, analytic expressions have been obtained for the coordi- 
nates along the intersecting M-line segment in terms of s, the arc length 
along the segment (equation (7.1.2)) and for the coordinates on the surface 
of the element in terms of the parametric variables u and w (equation 
(8.3.6)). Designating points on the surface by the subscript s and points 
on the M-line segment by the subscript £, expressions have been obtained 
for x^(s), y £ (s), z (s), x s (u,w), y g (u,w), and z s (u,w). At the intersec- 
tion point 

xjs) - x s (u,w) = 0 

y z (s) -y s (u,w) = 0 (8.3.20) 

z £ (s) - z s (u,w) = 0 

Equation (8.3.20) represents a system of three simultaneous nonlinear equa- 
tions in the three unknowns (u,w,s). 

The solution of this system of nonlinear equations requires an iterative 
procedure which must start from some estimate of the solution. An estimate 
of the solution, in terms of the coordinates of the intersection point 
(x. , y . , z i ) was obtained during the searching procedure to determine which 
M-line segment intersects which element of the intersected component. An 
estimate of s^ , the arc length along the M-line segment to the intersection 
point, was also obtained. The coordinate data must be used to determine esti- 
mates of the values of the u and w variables at the intersection point 
(u.j , w.j ). In the present method, the two planar subelements are converted 
into a single planar element and the values of u and w in the planar 
element are then determined. The planar element is formed from the three 
corner points defining the subelement in which the approximation to the inter- 
section point is located and a fourth point obtained by rotating the other 
corner point of the original element about the line separating the two sub- 
elements until it lies in the plane of the other three points. 
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The parametric coefficients of a point in a planar element bounded by 
straight lines can be determined most easily be first noting that all coef- 
ficients of second- and third-order terms in equation (8.3.6) must vanish 
and then by finding the other four coefficients directly from the properties 
of the sides of the element. The result is that equation (8.3.6) reduces to 

x(u.w) 3 Xgg + ( x i o -X 0(P U + ( x 01 “ X 0(P W + ^ x ll _x 10 - X 01 + x 00^ uw 

(8.3.21) 

Equations for thfe y and z coordinates have similar form. Given the 
values of x,y and z at a point, equation (8.3.21) and the corresponding 
y and z equations constitute a system of three equations for the two 
unknowns (u and w). The nonlinear term uw is eliminated from each equa- 
tion and then one equation is added to each of the other two equations to 
result, finally, in a system of two independent linear equations for the two 
unknowns . 

Given the estimate of the variables (s^ , u. , w.. ) at the intersection 
point, Newton's method is used to solve the nonlinear system (8.3.20). The 
process generally takes only four or five iterations to converge to an error 
in the square-root of the sum of the squares of the variables of less than 
10” 5 for a typical case. When solved, s i can be used in equation (7.1.2) 
and the corresponding equations for y and z to obtain the coordinates of 
the intersection point. Alternatively, u^ and w^ could be used in (8.3.6) 
and the corresponding y and z equations. 

8.3.5 Test Cases for the Intersection Method 

In order to verify the accuracy of the present method, a number of test 
cases were run. These are not examples of realistic aircraft components, but 
they were chosen because their intersection curves can be analytically deter- 
mined. Figure 18 shows a relatively simple case, the intersection of two 
circular cylinders. Using only four elements to represent each quadrant of 
each cylinder produces remarkably good results. There is no discernable dif- 
ference between the theoretical and the calculated results (to the scale 
plotted). This is an especially accurate case for any method, since the 
intersecting M-lines are straight lines and the elements on the intersected 



component have curvature in only one direction. A more difficult case, the 
intersection of two spheres, is shown in figure 19. This case tests the full 
capability of the intersection method. The intersecting M-line segments are 
circular arcs. The intersected elements have curvature in both directions and 
some of the intersected elements have zero-length sides (i.e., some elements 
are triangular). Nevertheless, the present method's results, using 250 ele- 
ments on each sphere, are very accurate; the calculated intersection curve is 
only very slightly different from the theoretical one. Figure 20 shows the 
final test case, the intersection of two ellipsoids. Like the case of the 
two spheres, this case shows good results, using 250 elements per body. The 
method failed to find intersection points for the M-lines parallel to the 
x-axis in figure 20, since these M-lines meet the intersected component 
tangential ly. 

In order to illustrate the method for a configuration more typical of 
aircraft configurations, figure 21 is included. This shows a wing intersect- 
ing a cylindrical section similar to the midsection of a typical fuselage. To 
give the intersection curve more character, the thickness of the wing has been 
greatly exaggerated. 



Figure 18. Intersection method test case — intersection of two circular 
cylinders (side view). 
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THEORETICAL INTERSECTION CURVE 

O CALCULATED INTERSECTION POINTS 

Figure 19. Intersection method test case 


y 


Figure 20. Intersection method test case — intersection of two ellipsoids 
(side view). 


THEORETICAL INTERSECTION CURVE 

° CALCULATED INTERSECTION POINTS 

— intersection of two spheres. 
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9.0 FINAL REPANELING OF COMPONENTS 
9.1 General Considerations 

The final major operation performed by the geometry package is the 
redistribution of points on all components which intersect other components. 
The point spacing requirements of the potential -flow method for intersecting 
bodies are not as well understood as those for isolated bodies. The accuracy 
of surface-singularity type methods in regions near concave corners is a 
matter of dispute in two-dimensions (reference 10) and has not been exten- 
sively studied in three dimensions. Some obvious requirements can be identi- 
fied, however. One is that portions of components which fall inside other 
components should be eliminated, or in certain cases, designated as extra 
strips or ignored elements. Another is that the repaneling should not cause 
abrupt changes in the element spacing. Therefore, the entire component should 
be repaneled to produce a smooth transition from the paneling in the region 
of the intersection curve to the paneling in distant regions, rather than just 
the region of the intersection curve being repaneled. It is also possible 
that future potential -flow methods (such as the method of reference 5) will 
require the matching of the corners of the elements of adjacent components to 
eliminate or at least reduce the size of any gaps between elements. 

It is possible to envision innumerable different types of intersections 
between bodies, each of which would have its own special paneling require- 
ments. It is not possible to develop a method sufficiently general to deal 
with them all. However, plausible paneling schemes can be developed for cer- 
tain frequently occurring configurations, such as wing-fuselage or wing-pylon 
cases, so that these cases can be handled routinely. The present method 
divides all components which are involved in intersections into three distinct 
general categories and provides separate means of dealing with each of them. 
All intersecting components (thosewith M-lines which pierce other components) 
are repaneled one way, with only minor variations (whether they are lifting 
or nonlifting). Nonlifting intersected components are repaneled differently, 
and lifting intersected components are repaneled in a different way still. 

The methods used for these three categories are described in the following 
sections. 
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9.2 Intersecting Components 


Typical examples of intersecting components include the wing in a wing- 
fuselage case, the pylon in a wing-pylon case, the winglet in a wing-winalet 
case, etc. Although all these examples are of lifting intersecting components, 
they can also be nonlifting (struts, for example). Lifting and nonlifting 
intersecting components are repaneled in essentially the same manner. There 
is no redistribution of points along N-lines, only a redistribution of N-lines. 
In every case, a new N-line along the intersection curve is added. If the 
component is of lifting type and if the intersecting end of the component is 
designated as having an extra strip, then the entire area inside the inter- 
sected component (from the intersection curve to the end N-line) is made into 
a single strip by eliminating all intervening N-lines. If the component is of 
non lifting type, or it does not have an extra strip at the intersecting end, 
then the portion inside the intersected component is simply eliminated. In 
the simplest cases, no other redistribution of points or N-lines is necessary. 
In most cases, however, the remaining N-lines on the exterior of the component 
must be moved in order to avoid irregularities in the widths of the strips. 

This is done by simply scaling the specified distances between N-lines by the 
fraction of the span of the component which lies outside the intersected 
component. To be more precise, the specified distances mentioned above are 
either the normalized distances between the planes of the N-lines (for the 
planar-section mode of operation) or the fraction of the arc lengths along 
the M-lines (for the arc-length mode of operation). The fraction of the span 
of the component which lies outside the intersected component is the fraction 
of the arc length along the M-line under consideration (for the arc-length 
mode of operation) or the fraction of the arc length along the first M-line 
of the component (for the planar-section mode of operation). The numerical 
techniques required to move the N-lines in this manner involve only interpola- 
tion procedures or cubic curve-plane intersection procedures which were 
described in section 7 above. 

An example of the program capabilities for repaneling intersecting 
components is shown in figure 22. Figure 22(a) shows a wing-fuselage case 
with the elements required for the initial geometric representation. In 
figure 22(b), the isolated components have been repaneled, the intersection 
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curve has been made an N-line and all exterior N-lines have been redistributed 
to produce a smooth distribution. In addition, the area inside the intersected 
component has been made into an extra strip, but for clarity of presentation, 

the extra strip has not been shown in the figure. 

More complicated cases, in which a body intersects more than one other 

body, can also be treated by this geometry package. However, in these cases, 

the intersecting body must be divided into more than one component. For 
example, in the case of the wing-fuselage with tip- tank shown in figure 23, 
the wing is divided (at fifty-percent semi span) into two components, an 
inboard component which intersects the fuselage and an outboard component 
which intersects the tip-tank. Since the same algorithm was used for distri- 
buting points on the N-lines of each component, and since points on the N-lines 
of intersecting components do not get redistributed, there is no mismatch in 
the elements at the junction between the components. Because of the redistri- 
bution of the N-lines, however, there is a small difference in the width of 
the strips in the inboard and outboard portions of the wing. Cases in which 
a body completely pierces several other bodies can be handled by dividing the 
intersecting body into components in an obvious way. 



Figure 23. Final repaneling of intersecting components — wing-fuselage-tip-tank 
case. 
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9.3 Nonlifting Intersected Components 

Typical examples of nonlifting intersected components include the fuselage 
in a wing-fuselage case, the tank in a wing-tip-tank case, etc. The most common 
application is to the wing-fuselage case. Therefore, the repaneling method in 
this geometry package was conceived primarily with fuselages in mind, but it 
should also be useful for other types of cases. 

Three options are provided for the final repaneling of nonlifting inter- 
sected components. The simplest option is to do no repaneling at all, in which 
case the final element distributions are as shown in figures 22 and 23. If this 
option is used, a portion of the intersected component falls inside the inter- 
secting component. If the intersecting component is not very thick, an adequate 
potential -flow solution may still be obtained in spite of this. If the inter- 
secting component covers a significant area on the intersected component, then 
the elements inside the intersecting component should be designated ignored ele- 
ments (elements with no singularities and no boundary conditions). Since there 
is no repaneling, however, some elements are only partially covered and the deci- 
sion whether or not the elements should be ignored requires some user judgment. 

For this reason, no mechanism has been provided for automatically designating 
elements to be ignored elements on nonlifting intersected components. If the 
user desires to use ignored elements, he must repanel the configuration using the 
geometry package, punch the resulting coordinates on cards, and terminate execu- 
tion of the program. Then he must execute the program a second time, using the 
punched output from the preceding case as input data, designating the ignored 
elements himself with the appropriate flags, skipping the geometry package, and 
proceeding straight to the potential -flow calculations. 

With the first option, the element distribution on the intersected component 
is not influenced by the element distribution on the intersecting component. 
Experience with the method of references 1 and 2 indicates that there is no strong 
need to repanel fuselages is most wing-fuselage cases; accuracy of the potential - 
flow solution is usually adequate using a distribution of elements appropriate 
for an isolated fuselage case. However, in some cases even better accuracy is 
desired. Also the newer potential -flow methods, using surface doublet distribu- 
tions, may require the elimination of the small gaps which result when the edges 
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of adjacent elements are not carefully aligned. The second and third options 
(both very similar) repanel nonlifting intersected components to produce element 
distributions which depend strongly on the element distributions on the inter- 
secting components. 

The calculations in both options start by defining a leading-edge and a 
trailing-edge point on the intersection curve. The trailing-edge point is 
assumed to be the first or last point on the intersection curve, since this curve 
is an N-line of the intersecting component and points on N-lines are input start- 
ing at the trailing edge, working around the perimeter, and ending back at the 
trailing edge. The leading-edge point is found by searching for the point on the 
intersection curve which has the smallest projected distance in the axial direc- 
tion from the front of the interesected component. Because of this, nonlifting 
intersected components which are repaneled with one of these two options must be 
input starting at the upstream end, must have roughly streamwise M-lines, and 
must be paneled using the planar-section mode of operation. The axial direction 
is defined by the direction cosines of the vector perpendicular to the planes of 
the N-lines. Planar cuts are then made through the leading and trailing-edge 
points on the intersection curve, perpendicular to the axis of the intersected 
component. The points of intersection of these planes and the M-lines of the 
intersected component define new N-lines. All N-lines forward of the leading- 
edge point and all those aft of the trailing-edge point are then redistributed 
in a manner similar to the redistribution of the N-lines on intersecting compo- 
nents. Forward of the leading-edge point, this is done by scaling the projected 
distances in the axial direction from the front of the component by a factor 
equal to the ratio of the distance from the front of the component to the leading- 
edge point and the distance from the front of the component to the first N-line 
aft of the leading-edge point. Aft of the trailing-edge point, the redistribu- 
tion is done the same way, except that the back of the component, the distance 
to the trailing-edge point and the distance to the first N-line in front of the 
trailing-edge point are used. 

The two options differ only in the distribution of N-lines in the region 
between the leading and trailing-edge points on the intersection curve. One 
of the two options passes planar sections through each point on the inter- 
section curve and uses the points of intersection of these planes and the 
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M-lines of the intersected component to define new N-lines. This option can be 
used when an element distribution with no gaps between adjacent elements is 
desired. The other option passes planar cuts through every second point 
on the intersection curve, starting with the leading-edge point and working 
aft on both upper and lower surfaces of the intersection curve. Since each 
surface may contain either an even or odd number of points, it is possible 
that planes will be passed through two consecutive points at the back of the 
intersection curve. This option can be used when gaps between adjacent 
elements can be tolerated, but it is desired to keep their size fairly small. 

Having redistributed the N-lines, it is next necessary to redistribute 
the points on each N-line. The method for doing this is the same for both 
available options. First, a search of the N-line passing through the leading- 
edge point on the intersection curve is conducted to find the point closest 
(in arc length along the N-line) to the leading-edge point. This point is 
moved to coincide with the leading-edge point. The other points on this 
N-line are redistributed in a manner which makes the resulting arc-length 
distribution reasonably smooth and similar to the previous distribution. For 
points on the N-line before the point closest to the leading-edge point on 
the intersection curve, this is done by scaling the arc length along the 
N-line by a factor equal to the ratio of the arc length to the leading-edge 

point and the arc length to the point closest to the leading-edge point. For 

points on the N-line after the point closest to the leading-edge point, arc 
lengths are scaled the same way, but using arc lengths from the other end of 
the N-line. Having determined the arc-length distribution of the redistributed 
points, the coordinates are found by the same interpolation procedure described 
in section 7.1. Points are redistributed on N-lines forward of the leading- 
edge point using the same arc-length distribution as on the N-line which 
passes through this point. Between leading-edge and trailing-edge points, 
points on N-lines are redistributed in a similar manner by scaling the arc 
lengths by a factor determined from the ratio of the distances to the inter- 
section point and to the point on the same M-line as the point closest to the 

leading-edge point. Aft of the trailing-edge point, points on all N-lines 

are redistributed using the same arc-length distribution as the redistributed 
points on the N-line passing through the trailing-edge point. Because the 
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N-lines through the intersection points between the leading- and trailing-edge 
points are broken by the intersection curve, it is necessary to break the 
component into more than one component. Presently, the program breaks the 
component into four smaller components, one forward of the leading-edge point, 
one below the intersection curve, one above the intersection curve, and one 
aft of the intersection curve. 

A typical example of the results of this repaneling method is shown in 
figure 24. This case shows a wing-fuselage with the wing at zero incidence, 
located slightly above the axis of the fuselage. In this case, use has been 
made of the option which passes an N-line through every point on the inter- 
section curve. Figure 25 shows a similar case, using the same option. In 
this case, however, the wing has ten degrees of incidence, illustrating the 
point that, for such cases, the present method bunches points below the 
intersection curve and spreads them out above the intersection curve. In 
extreme cases, this may result in an unacceptable element distribution with 
many elements bunched tightly near the bottom of the component and very few 
elements near the top. However, most cases will probably be much less 
extreme than the one shown, for which the element distribution is adequate. 
Figure 26 shows a wing-fuselage case similar to the one of figure 24, but 
with the wing located at the same level as the axis of the fuselage. In this 
case, the option which passes N-lines through every second point on the 
intersection curve has been used. 

9.4 Lifting Intersected Components 

The most common example of a lifting intersected component is the wing 
in a wing-pylon or wing-with-strut case. 

As in the case of nonlifting intersected components, lifting components 
are treated three different ways. Again, the simplest option is to do no 
repaneling at all. In this case, part of the intersecting component covers a 
portion of the surface of the intersected component and the remarks of sec- 
tion 9.3 concerning accuracy of the solutions, the desirability of using 
ignored elements, and the mechanism for generating them again apply. 
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SIDE VIEW (fuselage only) 



Figure 24. Final repaneling of nonlifting intersected components - wing-fuselage 
case (zero-incidence wing, N-lines through every point of intersec- 
tion curve). 


SIDE VIEW (fuselage only) 



Figure 25. Final repaneling of nonlifting intersected components - wing-fuselage 
case (wing with 10° incidence, N-lines through every point on inter- 
section curve). 
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Figure 26. Final repaneling of nonlifting intersection components - wing-fuselag; 
case (zero-incidence wing, N-lines through every other point on 
intersection curve). 

In some cases the intersecting component covers a significant portion of the 
surface of the intersected component, indicating a need for the use of ignored 
elements, but the paneling of the intersected component is so coarse that large 
gaps would be created by their use. The geometry package provides an option whic! 
greatly reduces the size of these gaps. In this option a search is conducted to 
find the points on the intersection curve having the minimum and maximum projecte< 
distance in the axial direction on the intersected component from the first point 
on the first N-line of the component. The intersected component must be paneled 
using the planar-section mode of operation for this calculation to be performed. 
The axial direction is defined by the direction cosines of the vector perpendic- 
ular to the planes of the N-lines. New N-lines passing through the minimum and 
maximum points on the intersection curve and lying in planes perpendicular to the 


axial direction are then created. The N-lines on either side of the intersection 
curve are redistributed in a smooth manner using the same procedure as the redis- 
tribution of N-lines upstream and downstream of the intersection curve on a non- 
lifting intersected component. No redistribution of points on the N-lines is 
performed. The elements covered by the intersecting component are automatically 
designated as ignored elements by the geometry package and it is possible to 
proceed directly to the potential -flow solution without checking or changing the 
results of the geometry package. Figure 27 shows the results of this repaneling 
option for the wing in a wing-pylon case. 

When more accuracy in the region of the intersection curve is desired, 
a third option can be used, resulting in no gaps between elements. This 
option first repanels the intersected component in a manner nearly identical 
to the one described above. The only difference is that the N-lines immedi- 
ately on either side of the intersection curve do not actually pass through 
any of the points on the intersection curve. Instead, each is offset a 
distance equal to the thickness of the intersecting component (projected in 
the axial direction). The elements covered by the intersecting component are 
again designated to be ignored elements. Now the method adds a new nonlifting 
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Figure 27. Final repaneling of lifting intersected components — wing-pylon case. 
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component with elements designed to fill in the gaps around the intersection 
curve. The new component has only one strip (two N-lines). One of the 
N-lines coincides with the intersection curve. The other N-line follows the 
boundaries of the quadrilateral region containing the ignored elements. 

Figure 28 shows the resulting element distribution for a case with very sparse 
point spacing on the pylon. Because of the correspondence of the numbers of 
points on the pylon and on the boundaries of the region containing the ignored 
elements, the paneling produced only quadrilateral elements (except at the lead- 
ing and trailing edges of the pylon). If the numbers do not correspond as in 
figure 28, then one or the other of the N-lines on the new component contains 
repeated points and triangular elements are produced. Points are repeated in 
such a way that the length of each M-line on the new component is minimized. 
Figure 29 shows a case having more points on the pylon than on the boundaries of 
the region containing the ignored elements, showing the repeated points on one 
N-line and the resulting triangular elements. 

N-lines which bound extra strips falling inside another component are not 
moved during the repaneling of lifting intersected components. Because of this 
and because none of the options described above in this section changes the 


BOTTOM VIEW 



Figure 28. Final repaneling of lifting intersected components — wing-pylon case 
(sparse element distribution on pylon). 
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BOTTOM VIEW 



Figure 29. Final repaneling of lifting intersected components — wing-pylon case 
(more points on pylon than surrounding region on wing). 

point distributions on N-lines, the repaneling of a lifting intersected compon- 
ent does nothing to destroy any match-up of the elements along the curve of inter- 
section of the component with another component which it pierces. Therefore, one 
component can play the role of both an intersected and an intersecting component 
(although one component cannot intersect or be intersected by more than one 
other component). Figure 30 shows the final element distribution for a wing- 
fuselage-pylon case. In this case, the wing and pylon were first repaneled as 
intersecting components; then the fuselage and the wing were repaneled as inter- 
sected components. 
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Figure 30. Final element distribution on a wing-fuselage-pylon case. 
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10.0 CONCLUSIONS 


The geometry package described above provides a means for significantly 
reducing the effort required to prepare the input data for three-dimensional 
potential -flow calculations. Data may be input to the program using either 
the original input format of reference 4 or the format of reference 5. 

Geometric input data generally consists of the coordinates of sets of sparsely 
defined points. In most cases, the number of points input can be at least an 
order of magnitude less than the number of points required for the potential - 
flow calculations. Each component of the configuration is automatically 
paneled using one of several algorithms provided by the method. The number of 
algorithms provided and the generality of some of them provide the user with a 
great deal of flexibility in determining the character of the resulting element 
distributions. Curves of intersection of components are automatically calcu- 
lated and all intersecting components are repaneled in a suitable manner for 
the potential-flow calculations. In many cases the potential -flow calculations 
can be performed in the same computer run as the geometry package calcula- 
tions, without interruption. To allow for intermediate checking, however, and 
to provide for those cases which cannot be run completely without interruption, 
provision has been made for the punched output of the defining coordinate data 
at each stage of the geometry calculation in a format suitable for reinput to 
the program. 
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